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Sound  Propagation  In  the  Sea: 

ANALYTIC  VERIFICATION  OF  WAVEFRONT  SPREADING  FORMULATIONS, 
NON -RECIPROCITY  IN  RAY  REVERSAL,  AND  FNWF  COMPUTER  PROGRAM 


1.  Definitions  and  Formulations 

In  Figure  l we'aiepiet* a segment  of  wavefront,  ATB;-of  sound 
propagating  in  a vertical  plane  in  the  se^  having  originated  at  the  source 
S.  The  trace  of  a point  r on  the  wavefront  defines  a sound  ray,  St, 
whose  propagation  instantaneously  is  along  the  wavefront-normal  unit- 
vector  tt  which  makes  the  angle  ) with  the  horizontal.  This  arbitrary 
ray  left  the  source  with  angle  ? . 


Fig.  1 Definitions 
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The  unit  vector  in  is  wavefront  tangent,  clockwise  normal  to  tt. 

Linear  coordinates  along  the  ray,  Sr  , along  the  wavefront,  rB  , and  along 
the  straight  line  defined  by  the  direction  in  at  r , are  denoted  by  s,  m and 
n respectively. 

The  formula  for  integrating  the  trace  may  be  expressed  by  the 
equation 

cos  y 1 n\ 

tt  • v y - - — - — - - in  • vC  U) 

The  first  term  on  the  right  hand  side  is  a coordinate  curvature  term  in 
which  R is  the  local  radius  of  curvature  of  the  horizontal  coordinate.  The 
second  term  is  the  refraction  where  C is  the  non-uniform  speed  of  sound 
in  the  sea. 

Not  only  is  the  depth  of  the  sea  small  compared  to  the  radius  of 
the  earth,  but  it  is  also  small  compared  to  the  earth's  oblateness.  This 
is  why  R is  held  constant,  as  the  mean  radius  of  the  earth,  in  numerically 
tracing  rays  at  the  Fleet  Numerical  Weather  Facility.  Thus  curvature  of 
the  horizontal  coordinate  is  taken  into  account,  but  not  strictly:  it  is 
given  a constant  value. 

It  should  also  be  noted  that  the  straight  line  defined  by  m , and 
dashed  in  Fig.  1,  takes  on  curvature  when  the  sea  surface  is  depicted  as 
a straight  line. 

The  formula  for  integrating  the  wavefront  divergence  along  the 
trace  may  be  obtained*  by  combining  the  spreading  rate 

tt  • V L = L in  • V y (2) 

*"The  wavefront-divergence  factor  in  ray-intensity  integration",  M.  M.  Holl 
Meteorology  International  Inc.,  Project  M-140  Technical  Note  Two, 
Contract  No.  N62271-67-M-2000 , April  1967. 
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L L 
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with  the  refraction  rate 


It  • v y 


l 

C 


in  • v C 


(3) 


to  obtain 


tt  * V (tt  . V L ) = - ^ in  • V ( in  • V C)  (4) 

L denotes  the  specific  wavefront  length  per  unit  radian  of  emission  angle, 
y . Equation  (4)  may  also  be  expressed  in  scalar  form  by 


2 2 
51  L a c 


and  by 

52L  _ L 

as2  = “ c 

In  treating  the  horizontal  as  straight,  in  applying  Eq.  (6) , a 

2 2 

curvature,  relevant  to,  but  negligible  in,  the  term  d C/d n , is  induced 
into  the  straight  line  defined  by  in.  However  this  curvature  is  not  that 
which  is  consistent  with  the  treatment  adopted  by  Eq.  (1) . For  analytical 
purposes  consistency  requires  that  we  rederive  the  wavefront  divergence 
equation  by  substituting  Eq.  (1)  for  Eq.  (3)  and  subsequently  neglecting 
the  curvature  of  the  horizontal. 

The  operation  tt  • V on  Eq.  (2)  yields 


5 2C  _idC.lL 

^ 2 L d s d s 
dn 


tt  • v (tt  • v L) 


r 


;L)  (in  • v y ) r L (tt  . Via)  . V>  r L in  • v (tt  • v> ) - L (in  . v tt)  • v y (7) 


(tt  • v L)  ( 


Wo  substitute 


tt  • V in 


- tt  tt  • vy 


. v tt  m m • vy 


and  obtain 


tt  • V (tt  • V L) 


(tt  • v 


/l.)  (m  • vy)  - L (tt  • V y)“  + L ui  • v (tt  • vy ) - L ( in  . v>  ) 


/e  next  selectively  substitute  from  bus.  (1)  and  (d): 


tt  • V (tt  V L) 


1. Vf)2  - L j— 2 + £ in  • VC  j - Lit  • v 


i.  (tu- 


!2  - L«.  V I + i m • V c ! 


which  may  be  developed  into 


tt  • V (tt  • V L) 


j 

l , „v  t \ 2 cos  y m • VC_  cos  J_y 

m * V (m  • vt.)  - L j R G r2 


sin  y tt  • VL 
R L 


f 


This  version  is  consistent  with  tracing  the  ray  according  to  Eq.  (1).  If 
we  neglect  the  curvature  altogether  by  setting  R = oo  then  Eq.  (12) 
reduces  to  Eq.  (4)  just  as  Eq.  (1)  reduces  to  Eq.  (3), 


Integration  of  L along  a ray,  by  Eq.  (12)  or  (6),  from  some  point, 
requires  specification  of  L and  9 L/9  s in  order  to  proceed  from  that 
point.  In  beginning  at  the  source,  we  specify 

Ls  = 0 (13) 

which  is  consistent  with  representing  the  source  effectively  as  a point 
source.  We  also  specify 

(9L/9s)s  = 1 (14) 

for  normalization  of  L as  the  specific  wavefront  length  per  unit  radian  of 

emission  angle,  y . 

o 

In  the  next  section  we  demonstrate  by  analytic  example  that  the 
wavefront  divergence  generally  depends  not  only  on  the  path  of  the  ray 
between  source  and  receiver  but  also  on  which  end  of  this  path  is  the 
point  source.  Generally,  the  relative  intensity  loss  changes  in  inter- 
changing the  position  of  the  source  and  receiver. 

Another  way  of  looking  at  this  fact  is  by  the  statement  that  the 
proportional  change  in  intensity,  in  propagating  through  a finite  segment 
of  path  in  one  direction,  depends  on  the  curvature  of  the  wavefront  at 
the  commencement  of  the  path  segment.  The  integration  is  reversible  if 
we  duplicate  the  direct  exit  curvature  of  the  wavefront  for  initializing 
the  reverse  traverse  of  the  segment. 
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2. 


An  Analytic  Solution 


In  the  present  section  we  design  a parameterized  analytic 
solution,  showing  its  complete  derivation  from  Eq.  (1).  This  solution 
serves  to  (a)  test  any  of  our  derived  equations , (b)  edify  aspects  of  the 
propagation,  and  (c)  test  and  evaluate  computer  programs  of  the  Fleet 
Numerical  Weather  Facility  for  numerical  integration  of  Eqs.  (1)  and  (6). 

For  speed-of-sound  distributions  which  are  functions  of  depth, 

Z , only,  Eq.  (1)  reduces  to 


9V  _ _1  _ SK 
dx  R dZ 


(15) 


where  K = toC.  Our  solution  is  for  the  special  profile 


C = CQ  e "kZ  (16) 

and  a point  source  at  arbitrary  depth  Zg  . The  horizontal  (i.e.  range) 
coordinate  is  x _>  0 . We  shall  only  solve  for  direct  rays , not  going 
into  sea-surface  or  bottom  reflections.  The  constant  k may  be  positive 
or  negative;  we  consider  it  to  be  positive  for  discussion  purposes. 

The  Eq.  (16)  profile  submitted  to  Eq.  (15)  yields 


ay 

dx 


(17) 


where  /i  s k - l/R  is  considered  to  be  a positive  constant.  Integration  of 
Eq.  (17)  yields 


y = ys  + 


HX 


(18) 


-6 


for  any  and  all  rays.  All  rays  are  refracted  downwards,  with  maximum 
range 


x = - (tt/2  - yJ 
m n S 


(19) 


as  asymptote.  The  angle  at  which  a ray  leaves  the  point  source  is,  by 
definition,  limited  to  the  range 


- rr/2  < yg  < n/2 


for  a ray  defined  by  yg  we  have 


d2 

“ = tan  y = tan  (yg  + ^x) 


(20) 


(21) 


Integration  yields 


Z - Z„  = - - 


cos  (y^  + ^x) 
cos  y„ 


(22) 


for  any  and  all  rays.  This  gives  the  limiting  ray,  which  grazes  the 
surface  Z = 0 at  angle  T =0,  to  be 


= arc  cos  e 

b # 1 


- 


(23) 


The  maximum  range  for  direct  rays  is  thus 


1 ,v  - /iZ_ . 

x = — (t  - arc  cos  e n S) 
M n 2 


(24) 


i 


I 


% 


Equation  (22)  may  be  transformed  into 


un>s 


- ii  (Z  - Zg) 

cos  fix  - e 

sin  |ix 


(25) 


which  gives  y at  any  point  within  direct  range.  It  is  pertinent  to  note  that 
the  distribution  is  single  valued.  This  means  that  no  rays  cross.  There 
are  no  caustics  — a caustic  being  the  locus  of  a wavefront  folding  point. 

By  definition 


V y 2 ui  I (26) 

where  in  is  defined  in  Fig,  1 and  L is  referred  to  as  the  specific  wave- 
front  length  per  unit  radian  of  emission  angle.  The  ascendent  operation  on 
Eq.  (25)  yields 


-fi(Z"Z  J 


sec  >s  vyg  = 


cos  ux  e 


|i  vx 


sin  /ix 


- M (Z-ZQ) 


sin  fix 


v Z 


(27) 


where  vx  and  VZ  are  orthogonal  unit  vectors.  The  magnitude,  with  sign 
consideration,  leads  to 


-M(Z-Zg)  -2ji  (Z-Zg) 

1-2  cos  fjx  e + e 


1/2 


(28) 
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This  is  our  analytic  solution. 


Equation  (28)  may  also  be  transformed  into 


L 


1 1/2  2 (Z  ZS) 

- V e 

M 


(29) 


where 


p e 2 cosh  p (Z  - Z ) - 2 cos  px 


(30) 


and  which  shows  that  L is  positive  everywhere. 

For  later  usefulness  we  combine  Eqs.  (25)  and  (28)  to  form 


L ■ 


sin  ux 
H cos  >g 


(31) 


We  also  note  that 


tan  y = 


fi(Z  - Zg) 

e - cos  ux 

sin  ux 


(32) 


For  an  arbitrary  wavefront  in  the  x , Z plane , 


dx  . . _ 

- dZ  ■ tan  y - 


|l(Z  - z ) 

e - cos  ux 

sin  fxx 


(33) 


which  integrates  to  yield  the  wavefront  curve 
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m(Z-Zs)  - p (Z  - Zg ) 

2 cos  px  * e " + A c 


The  value  of  A determines  the  particular  wavefront:  for  the  wavefront 


passing  through  the  point  x*  , ZR  , 


A - 2 cos  pxA  - 1 


Non-reciprocity  by  Analytical  Solution 


For  an  effective  point  source  at  arbitrary  location,  x 0 and 


depth  Z.  , Eg.  (29)  gives  the  specific  wavefront  length,  L,  at  an 

u 


arbitrary  receiver  location  within  range,  at  xR  and  depth  ZR  : 


\ 1/2  “ *2  (ZR  “ ZS ) 

LR  “ M ® 


where 


Od  e 2 cosh  n (ZR  - Zg)  - 2 cos  px 


The  particular  direct  ray,  making  the  connection,  is  given  by 


tan  V 


cos  |ixD  - e 


-"(ZR-ZS) 


sin  pxR 


If  we  now  interchange  locations  of  source  and  receiver  then  the  specific 


wavefront  length  at  Z , by  ray  traversing  the  same  path  as  before  but 

o 
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in  the  reverse  direction,  is 


1 1/2  " 2 (ZS  " V 

Ls  - ; °s  e 


where 


Og  = 2 cosh  4 (Zg  - ZR)-  2 cos/ixR 


= Or 


(39) 


(40) 


We  note  that  <p  is  an  even  function  but  the  exponential  factor  in  L is 
not.  The  ratio  of  the  spreading  amounts  is 


(41) 


In  our  example  the  downward  traverse  of  a path  produces  less  spreading 
than  the  upward  traverse  of  the  same  path. 

We  have  shown  by  analytical  example  that  the  field  of  listening 
coverage  for  a receiver  location  cannot  be  obtained  on  the  basis  of 
equivalency  by  treating  the  receiver  as  a normalized  source.  However 
while  not  strictly  equivalent  such  procedure  may  yield  reasonable 
approximations  under  most  conditions. 

We  can  edify  the  non  reciprocity  in  a more  general  way.  For 
all  ray  paths  along  which 

-^7  = 0 (42) 

dn 
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Eq.  (6)  reduces  to 


9L  1 ac 
9s  C 9 s 


Written  as  perfect  differentials,  this  yields 


where  Cg  is  the  sound  speed  at  the  source  location.  In  integrating  along 
such  a ray  path  from  a source  at  location  A and  a receiver  at  location  B, 


we  obtain 


C ds 


Along  the  same  path,  for  a source  at  location  B and  a receiver  at  location 
A,  we  obtain 


C ds 


The  ratio  for  alternate  directions  of  traverse,  yields 


is  £b 

LA  " CA 


In  this  regard,  intensity  is  maximized  toi  a lecuivei  located  at  a sound- 
speed  minimum. 


-1.  Analytical  Test  of  Spreading  Equation 

Out  analytic  solution  was  completely  derived  Horn  Eq.  (1)  with 
R rt  constant.  This  solution  may  be  used  as  a test  of  other  equations  of 
Section  l aiui  any  further  derivations.  We  will  use  it  to  test  l.q.  (12) 
which  reduces  to  Eq.  (t»)  for  R - We  express  Eq.  (12)  in  the  form 


2 

a"i. 

9- 


L i ifc  1 9C  Jkj 

0 \ bnZ  - l c«s  | 


L 


I 2 cos  ^ a C*  cos"  y sin  y Al.  ^ 
j R C an  1 r2  ~ RL  as  | 


(4tl) 


The  test  we  make  is  with  L in  the  form  expressed  by  Eq.  (31).  Noting 
that 


a s 


o 


(-10) 


ax 

as  ^ as 


p cos  y 


(so) 


we  obtain 


al  = iva  ms&z-l  _ ui  ^2  (sl) 

cAs  cos  yc  H tan  jix 

o 


(52) 


as 


2 


- L 


2 

cos 


y 


2 sin  y cos  y 

^ tan  fix 


» - k C cos  y (53) 

a n 


. 2 


+ k“  C cos  2 > 


(54) 


* - k C sin  y (55) 

o s 

Substitution  of  Eqs.  (51)  through  (55)  In  Eq.  (48),  with  p = k - lA, 
readily  shows  identity. 


5.  A Test  of  the  Computer  Programs 

In  applying  the  analytic  solution  of  Section  2 as  a test  for  the 
computer  programs  of  the  Fleet  Numerical  Weathei  1 aclllty  we  must 
specify  numerical  values.  For  the  sound-speed  profile,  Eq.  (lb),  we 
specify 

. -1 

c - 1,680  meters  seconds 

-5  - 1 

k * 2 x 10  meters 

and  submit  the  sound-speed  profile  values: 


V 
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z (meters)  C (meters  seconds  1 ) 


0 

1,680.00 

500 

1,663.28 

1000 

1,646.73 

1500 

1,630.35 

2000 

1,614.13 

2500 

1,598.06 

3000 

1,582.17 

3500 

1,566.42 

•1000 

1,550.83 

•1500 

1,535.40 

5000 

1,520.  13 

5500 

1,505.00 

6000 

1,490.03 

6500 

1,475.20 

7000 

1,460.52 

7500 

1,445.99 

8000 

1,431.60 

8500 

1,417.36 

Tests  one  and  two  were  run  without  the  curvature  factor.  That 
is,  R = , and  the  rays  are  obtained  by  integrating  Eq.  (3)  and  the 

spreading  factor  by  Eq.  (6). 


In  Test  one  the  source  is  specified  at  the  1000  meter  depth.  The 
bottom  is  specified  at  8,500  meters.  The  rays,  T , from  - 10°  to  + 10'' 
in  two-degree  intervals,  are  requested.  The  integration  of  each  ray  is 
terminated  at  the  first  bottom  intercept.  From  the  output  we  select  the 
results  shown  in  Table  I;  we  also  show  evaluations  of  Eqs.  (25)  and  (28) 
for  the  same  locations. 
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Table  I:  No  Curvature 


I'.qs.  (3)  and  (t>) 

Integrated 

Eqs.  (25)  and  (28) 
Evaluated 

X (yds) 

Z (feet) 

V 

S 

(degrees) 

L (yds) 

>^>  (degrees) 

L (yds) 

40, 109 

27,887.14 

- 

10.00 

37,220.0 

- 9.99 

57,l/0 

35,437 

27,887.14 

- 

6.00 

33,175.4 

- 6 . 00 

33,190 

31,  174 

27,887. 14 

- 

2.00 

29,528.9 

- 2.00 

29 , 53  j 

29 , 209 

27,887.14 

0 . 00 

27,841.3 

0.00 

27,840 

23,984 

27,887. 14 

6.00 

23,352.4 

6.00 

23-356 

21,039 

27,887.14 

10.00 

20,841.3 

10.00 

20 , 84q 

1,909 

3,180.87 

- 

2 . 00 

1 , 909 . 7 

- 2.00 

1,910 

Test  two  has  the  source  at  the  bottom,  8,500  meters,  and  selects 
rays,  > , as  the  reflected  bottom  intercept  angles  of  Test  one.  We  show 

only  a result  tor  the  reversal  ot  the  zero  ray  of  Test  one  for  which  we 
evaluate  1.  by  Tin.  (31)  tor  the  given  V : 

s (yds)  Z (feet)  >s  (degrees)  L (yds)  Tig.  (31) 

29,209  3,280.84  - 30.60  32,362.9  32,344 

In  the  opposite  propagation  the  spreading  for  this  ray  path,  given  by 
Test  one,  amounted  to  27,841  yards. 

Test  three  repeats  Test  one  but  with  curvature;  K (>,371  km, 

— ^ 

making  p 1.9843  x 10  . The  numerical  program  obtains  the  rays  by 

integration  of  i'.g.  (1)  and  the  spreading  by  integration  of  Kq.  (6) . Our 
analytical  solution,  however,  satisfies  Eq.  (12).  This  inconsistency 
in  the  treatment  of  curvature  is  included  with  truncation  errors  in  the 
disparity  between  the  integrated  results  and  the  evaluation  of  the 
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analytic  solutions.  The  selected  results  are  shown  in  Table  II. 


Table  II:  With  Curvature 


Eqs.  (1)  and  (6)  Integrated 

Eqs.  (25)  then  (31) 
Evaluated 

x (yds) 

Z (feet) 

>g  (degrees) 

L (yds) 

(degrees) 

L (yds) 

40,323 

27,887.14 

- 10.00 

37,401.2 

- 9.99 

37 , 387 

35,609 

27,887.14 

- 6.00 

33,304.0 

- 6.00 

33,365 

31,311 

27,887.14 

- 2.00 

29,633.7 

- 2.00 

29,672 

29,330 

27,887.14 

0.00 

27,934.9 

0.00 

27,965 

24,066 

27,887.14 

6.00 

23,416.7 

6.00 

23,436 

21,102 

27,887.14 

10.00 

20,890.2 

10.00 

20 , 90q 

1,924 

3,180.08 

- 2.00 

1,924.8 

- 2.00 

1 » 9 2 5 

These  results  check  out  the  particular  FNWF  computer  programs  and 


attest  to  their  accuracy. 


6.  A Note  on  the  Lateral  Spreading 

In  general  a wavefront  is  a surface  having  two  dimensions  of 
spreading.  The  specific  wavefront  area  may  be  expressed  by 

E = H L (56) 

where  H is  the  specific  wavefront  length  measured  horizontally  and  L is 
the  specific  wavefront  length  measured  normal  to  H.  Both  H and  L may 
be  integrated  along  a ray  according  to 


L, 
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. A A 


j 


(57) 


a2H 

H I 

1 afc 

I 

ac 

aH 

*s2 

c I 

1 ah2 

H 

a s 

a s 

a2L 

L J 

i a2c 

\ 

ac 

aL 

^ ' 

C 1 

t an2  " 

L 

as 

a s 

(58) 


where  h is  a linear  coordinate  along  the  horizontal  normal  to  the  ray  and 
n is  normal  to  coordinates  h and  s. 

Sections  1 through  5 deal  with  ray  tracing  in  a vertical  x,Z  plane 
and  with  the  L component  of  the  spreading.  That  the  rays  stay  in  the 
plane  implies  the  condition,  or  approximation,  that 


ac 
a y 


H 


0 


(59) 


where  y is  normal  to  the  x,Z  plane.  If  we  further  assume  that 

a2c  afc  9 

ah2  ay2 


(60) 


then  Eq.  (57)  reduces  to 


a2H  aH  i ac 

.2  a s c a s 

a s 


(61) 


Eq.  (61)  is  a perfect-differential  relationship  which  yields 


aH  = _C_ 

as  cs 


(62) 
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where  C is  the  speed  of  sound  at  the  source  location.  The  integration 
b 

of  H takes  the  form 


(63) 


Again  we  note  that  there  is  no  reciprocity  in  interchanging  source  and 
receiver  locations  for  a specific  ray  path.  The  ratio  of  spreading  amounts 
in  the  two  directions  is  inverse  to  the  ratio  of  source  speeds , again 
implying  best  listening  at  speed-of-sound  minima. 

For  speed-of-sound  as  a function  of  depth  only,  we  have  Snell's 

law: 


cos  y 
cos  Tg 

Substitution  in  Eq.  (63)  yields 

H = -r- 

COS  7g 

A generalization  of  the  last  expression  may  be  obtained  for  radial 
symmetry  in  the  speed-of-sound  distribution;  the  symmetry  is  that 
expressed  by  x 2 R,  the  horizontal  range  in  all  directions  from  the  source 
location. 


_C_ 

C„ 


(64) 


^ cos  y ds 
0 


cos  X, 


(65) 


r 
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